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Abstract 

Time series of heartbeat activity of humans can exhibit long-range correlations. 
In this paper we show that such kind of correlations can exist for the heartbeat 
activity of much simpler species like Drosophila melanogaster. By means of the 
method of multifractal detrended fluctuation analysis (MFDFA) we calculate frac- 
tal spectra /(a) and h{q) and investigate the correlation properties of heartbeat 
activity of Drosophila with genetic hearth defects for three consequent generations 
of species. We observe that opposite to the case of humans the time series of the 
heartbeat activity of healtly Drosophila do not have scaling properties. Time series 
from flies with genetic defects can be long-range correllated and can have multifrac- 
tal properties. The fractal heartbeat dynamics of Drosophila is transferred from 
generation to generation. 

1 Introduction 

The irregular and complex structure of the time series (ECG) of human heartbeat dynam- 
ics is an object of considerable clinical and research interest pQ, This structure is 
connected not only to the external and internal perturbations but also depends on the syn- 
ergetic action of muscle and nervous systems which influences the correlation properties 
of the time series. In many simple systems the correlation function of the measured time 
series usually decays exponentially with the time. In complex systems the correlations 
can decay with power law and because no characteristic scale is associated with the power 
law such systems are called scale-free. Their correlations are called long-range because at 
large time scales the power law function is always larger than the exponential function. 
Below we are interested in the presence of long-range correlations in the time series for 
heartbeat activity of Drosophila melanogaster - the classical object of Genetics. Due to 
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the short reproduction cycle of Drosophila we can investigate the correlation properties 
of the heartbeat dynamics for three consequent generations. This allows us to study the 
relation between genetic properties of Drosophila and correlation properties of the time 
series of its heartbeat activity. 

The paper is organized as follows. In Sect. 2 we describe the investigated system, 
recording of the time series and quantities used for their analysis. The analysis of the 
obtained fractal spectra is performed in Sect. 3 . Some concluding remarcs are summarized 
in the last section. 

2 System and methods 

We investigate time series of the heart activity (ECG) of Drosophila melanogaster ob- 
tained from mutant flies and wild type controls provided by Bloomington Drosophila 
Stock Center, U.S.A. We crossed male Dopa decarboxilase (Ddc) mutant (FBgn 0000422 
located in chromosome 2, locus 37C1) and female shibire (shi) (FBgn 0003392 located in 
chromosome 1, locus 13F7-12). The Ddc mutants' heartbeat rate is about 60 % of the 
normal one. Ddc codes for an enzyme necessary for the synthesis of four neurotransmit- 
ters: norepinephrine, dopamine, octopamine, serotonin, related to learning and memory. 
The shibire (shi) mutants cause paralysis at high temperature. They code for the pro- 
tein dynamin, necessary for the endocytosis. Its damaging at high temperature stops the 
transmission of the impulse through the synapses, causes paralysis, and eliminates the 
effect of the neurotransmitters on the heart |3]. ECGs were taken from three consequent 
generations of species. Drosophila heartbeat was recorded optically and digitalized. Opti- 
cal ECG records were taken at a stage PI (white puparium) of a Drosophila development 
when it is both immobile and transparent and the dorsal vessel is easily viewed. The 
object was placed on a glass slide in a drop of distilled water under a microscope (magni- 
fication 350 x). Fluctuation in light intensity due to movement of the dorsal vessel tissue 
was captured by photocells fitted to the one eyepiece of the microscope. The captured 
analogue signal was then digitized at 1 kHz sampling rate by data acquisition card and 
Lab VIEW data capturing software supplied by National Instruments. 600000 data points 
were taken for each sample. 

The obtained time series are analysed by the multifractal formalism which is widely 
used in mathematics, physics, and biology [T], [S], [H], [7], jH]. The investigation is based 
on the spectrum h{q) of the local Hurst exponent and on the fractal spectrum /(a) piTUj. 
Let us consider a set of points which lies in an A^— dimensional Cartesian space covered 
by a grid of A^- dimensional cubes of edge length e. If for small e we need A^*(e) cubes to 
cover our set we can define the spectrum of generalized dimensions 



where g is a continuous index, /ifc is the natural measure, i.e., it is a measure of the fre- 
quency with which a typical orbit visits various cubes covering the investigated attracting 
set of points for the limit case when the length of the orbit goes to infinity (in addition 
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the frequences have to be the same for all initial conditions in the basin of attraction of 
the attractor except for a set with Lebesque measure 0). Thus for /i^ we have 

^(cfc,xo,T) 

lik = hm — , 2 

where ^ is the time the orbit originating from xq spends in the cube Ck in the time interval 
< t < T. -Do is called capacity of the set and it is not integer for some sets. From (^J 
by means of the L'Hospital rule we can easily obtain 

D, = lim^£i^^l/^ (3) 
e^o Ine 

Di is called also information dimension (as it measures how the information is scaled with 
ln(l/e)). In general > Di > D2 > . . .. If Dq varies with q the measure, associated 
with Dq is called multifractal measure. 

Let a set S be covered with a grid of cubes of unit size e and /i is the probability 
measure on S (/u(S') = 1). Let /i(cjt) = /ifc where Ck denotes again the k—th cube. We can 
assign a singularity measure to each cube 

fik = e"'^ (4) 

For small e we can make continuous approximation for the number of cubes for which 
is between a and a + da, i.e., we can denote this number as p{a)e~-^^°'^da. Substituting 
(jlj in the relationship for I{q, e) and after a transition from a sum over the cubes to an 
integration over the a we obtain 

N*{e) 



I{q^e) = X; e"**? = / da*p{a*)e~f^''*h'^''* = 
k=i 

= J da*p{a*) exp {[/(a*) - qa*] ln(l/e)} (5) 

For small e ln(l/e) is large and the main contribution to the above integral is from the 
neighborhood of the maximum value of the f{a*) — qa*. Let f{a*) be smooth. The 
maximum is located at a* = a{q) given by 

^ :[/(«*) - \a*=aig)= ^ \a*=a= q (6) 



da*'-" ' -.-J, "-v.; 

[/(a*) - \a*=a{q)= ^ \a*=a= q (7) 



Now we take the Taylor series representation of the function F{a*,q) = f{a*) — qa 
around the point a* = a{q) and substitute it in The result is 

I{q,e) = exp{[f{a{q)) - qa] ln(l/e)} x 
X J rfa*p(a*)e-(i/2)/"("W)(-*-"W)' 

^exp{[/(a(g))-ga]ln(l/e)} (8 
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and a substitution of relationship (jH)) in ((H) leads to 

1 



I - (9) 



Using (jni) we obtain 



Then 



^[(.-1)Z.,]^«(.)^| (10) 
T{q) = {q-l)D,^D, = ^ (11) 



From © 



/(«(?)) = - - 1)^, = - ^(^) (12) 

For each g from (fTT|) and ()12p we can obtain a{q) and /(a) thus parametrically specifying 
the function /(a). And a can be connected to the local Hurst exponent by means of the 
relationships 

a = h{q) + q^, f (a) = q[a - h{q)] + 1 (13) 
Thus obtaining the h{q) spectrum we can obtain also a and /(a) spectra by means of 

For calculation of h from the heartbeat time series we can use the method of multifrac- 
tal detrended fluctuation analysis (MFDFA) or the more complex wavelet transform mod- 
ulus maxima method (WTMM), initially developed for investigation of quasi-singularities 
of turbulent signals (for applications of this method see [HI [12], [IH|, [Hj, [IH], [il^). In 
this paper we shall use the MFDFA method which realization is as follows ^Zj. First of 
all we have to calculate the profile function Yi. For this we calculate the mean (x) of the 
investigated time series {xk} and use it to obtain the profile function 

Y, = j2ixk-{x)), t = l,2,...,N. (14) 
fc=i 

The following step is to divide the time series into segments and to calculate the variation 
for each segment. The division is into Ns =int(A^/s) segments and because the obtained 
segments would not include some data at the end of the investigated time series, additional 
Ns segments are added, which start from the last value of the sequence in the direction 
to the first value of sequence. 

In order to calculate the variation we have to calculate the local trend (the fitting 
polynomial yu{i) for each segment of length s, where s is between an appropriate minimum 
and maximum value). Then the variations are defined as 

F'{u, s) = - f: {Y[iu - l)s + ^]~ y,i^)}' (15) 
for the first Ng segments and 



F\i^, s) = -j2{y[N-{i^- N)s + - y,{z)}' (16) 
^ i=i 
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for the second Ng segments. Finally we construct the g-th order fluctuation function 

F,{s) = {[l/{2Ng)]Y.[F\iy,s)]^/'y/''. (17) 

The scaling properties of Fq{s) determine the kind of fractal characteristics of the time 
series. For monofractal time series Fq{s) scales as s of constant power h for each q. For 
sequences of random numbers this constant h has the value 1/2. Even in presence of 
local correlations extending up to a characteristic range s* the exponent h = 1/2 would 
be unchanged when s » s*. If the correlations do not have characteristic lengths the 
exponent h would be different from 1/2. 

The procedure described above is appopriate for determination of positive Hurst ex- 
ponents which are not very close to zero. For close to zero or negative exponents we have 
to add a step after the calculation of the profile function namely to calculate the profile 
function of the profile function Y 

Y: = t[Yik)-{Y)] (18) 

k=l 

and the function Y* should be used further in the MFDFA procedure. The result is 
that if there is a scaling in the fluctuation function this scaling is connected to the Hurst 
exponent as 

F;(s) oc /('?)+^ (19) 

In our investigation below we use MFDFA(l) i.e. the local trend for each segment is 
approximated by a straight line. 

3 Results 

A part of typical time series for the heart activity of Drosophila melanogaster is presented 
in panel (a) of Fig. 1. From these time series we can construct time series for the 
interbeat intervals (presented in panel (b) of Fig. 1). Such time series are widely studied 
for humans ^H], ^H] because they can be easily measured in a noninvasive way and 
may have diagnostic and prognostic value. The interbeat time series of human heartbeat 
dynamics has (i) monofractal properties (constant h) for humans with heart diseases and 
(ii) multifractal properties (nonconstant h) for time series from healtly humans. As we 
shall see this is not the case for Drosophila. 

The autocorrelation functions for a healtly control fly and for the parents with heart 
defects are shown in Fig. 2. In all three panels we observe that a significant degree of 
correlation exists even for large values of n. In addition in panel (c) we observe systematic 
decrease of the autocorrelation function and transition from predominantly correlated 
behavior for small n to predominantly anticorrelated behavior for large n. Thus the 
dynamical consequences of the different genetic heart defects of Drosophila are clearly 
visible. 

Panel (a) of Fig. 3 shows the fiuctuation functions ( g = 2) for a healtly Drosophila 
and for parent files with heart defects. In the case of humans the normal sinus rhythm of 
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Figure 1: Panel (a): Typical time series of the heart activity of Drosophila melanogaster. 
The unit for time is 0.001 s. Panel (b): Interbeat intervals for the time series of heart 
activity of Drosophila. As we can see the time series of the heart rate fluctuate irregularly 
from beat to beat. 



the heartbeat activity has complex behavior similar to the behavior of a chaotic attractor 
[I]. The heart dynamics of humans with heart diseases may become more periodic in 
comparison to the heartbeat dynamics of the healtly individuals. The heartbeat dynamics 
of the investigated here Drosophila shows opposite behavior. We see that the fluctuation 
function for the healtly Drosophila does not exhibit scaling at least for small s and this 
lack of scaling is observed for all values of the parameter q. The deviation from the scaling 
behavior for the fluctuation function means that we can not calculate any fractal spectra 
for the healtly Drosophila opposite to the case of the flies with genetic defects where the 
fluctuation function can show good scaling properties for the whole studied range of s. 
We note that the fluctuation functions for the parents seem to be very close to a straight 
line on a log-log scale. Thus we shall proceed with calculation of the fractal spectra. 
These spectra will have different properties for time series of Drosophila with different 
heart defects. 

The Hurst exponent for the two parents is presented in Fig. 4. h is not a constant and 
hence the two time series of the parents have multifractal properties. Thus multifractal 
cardiac dynamics can be observed not only for humans but also for much simpler animals 
like Drososphila. In figures 5 and 6 we see the kinds of spectra of the Hurst exponent char- 
acteristic for the first and second generations of flies obtained from the above-mentioned 
parents with genetically defect hearts. The spectra in panels (a) and (b) in Fig. 5 are 
of the same kinds as the spectra of the two parents. The spectrum in panel (c) has non- 
typical from the point of view of physics because in most physical systems h decreases 
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Figure 2: Autocorrelation function C{n) for the time series of the heart activity of 
Drosophila. Panel (a) : autocorrelation for a healtly animal. Panels (b) and (c): au- 
tocorrelations for the two parents: female (panel (b)) and male (panel (c)). 



with increasing q. For the second generation of flies we observe the two kinds of h{q) 
spectra existing in the case of the parents plus an additional kind of spectrum with Hurst 
exponent which is systematically smaller than 0.5 for positive q i.e. the anticorrelations 
dominate the corresponding time series. 

The difference in the dynamical properties of the intermaxima time series for the 
heartbeat activity of Drosophila can be investigated by means of their f{a) spectra. Fig. 
7 shows these spectra for the parents. For the spectra with parabohc form, the parts 
of elements of the time series with a given value of a, build a partial fractal with a 
fractal dimension denoted by f{a). The top part of the spectrum which is located around 
some value /(a*) corresponds to the statistical most significant part of the spectrum 
(corresponding to the parts of the time series with the largest dimension), /(a*) gives 
the value of this largest dimension and we can distinguish the time series with respect 
to the value of a* and the width of the spectrum around the maximum ( A = ar{f*) — 
ai{f*), where /* is characteristic which we shall take to be equal of 0.9 f max in order to 
compare the parameters of the f{a) spectra of all generations of Drosophila. ai and 
are the values of a corresponding to /* and positioned to the left and to the right with 
respect to the value a* corresponding to the maximum of the f{a) spectrum). Wide f{a) 
spectrum corresponds to more distributed multifractal (the partial fractal dimensions are 
less concentrated around the maximum partial dimension fmax) and a narrow spectrum 
corresponds to more concentrated multifractal. Coming back to the spectra of parents in 
Fig. 7 we observe the typical parabolic form of the spectrum only for the male parent. 
Thus the form of the f{a) spectrum can help us to distinguish among the heart defects 
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Figure 3: Panel (a): Typical forms of the fluctuation function for intermaxima time series 
with and without scaling properties. Fluctuation function ^2(5) for the parents are as 
follows. Circles: female parent. Squares: male parent. For these time series r.m.s. fit of 
the power law is shown as a continuous straight line. The closeness to a straight line on 
the log-log scale means that the corresponding time series of the intermaxima intervals 
have scaling properties. For comparison typical time series for a healtly Drosophila (filled 
triangles) is presented. We do not observe scaling and thus we cannot calculate any fractal 
spectra. Panel (b): Fluctuation functions -^2(5) and power- law r.m.s. fits (sohd hues) for 
time series of the first generation of flies (the kids) . As we see there is no drastic breaking 
of the scaling as it is for the healtly Drosophila of panel (a). 

of Drosophila as some of these defects (and in particular the genetic defect of the female 
parent) can lead to nonparabolic form of the /(«) spectrum, i.e., to deviation from the 
ideal multifractal behaviour, fmax = 1 for the spectrum of the male parent and its 
0.9/maa: width is A = 0.18. The result of the combination of the two kinds of dynamics 
leading to parabolic and nonparabolic f{a) spectra can be observed in the spectra of the 
two generation of flies following the parents. The characteristic spectra for the second 
generation are presented in Fig. 8. We observe two kinds of consequences from the form 
of the spectrum of the female parent (i) the nonparabolic kind of spectrum is reproduced 
as it can be seen in panel (b) of Fig. 8. and (ii) some (but not al) of the parameters of 
the parabohc spectra change. We note that for the parabolic spectra in panels (a), (c), 
(d) of Fig. 8 fmax remains unchanged and equal to 1 not only for this generation of flies 
but also for the parabolic spectra in the next generation shown in Fig. 9. For the second 
generation of flies a for f^nax is dispersed around 0.66 - its value for the male parent. 

In the third generation of flies the nonparabolic form of the spectrum is reproduced 
again. From several characteristic examples of parabolic spectra of this generation which 
are shown in Fig. 9 only one of the spectra has a wide basis. For all spectra fmax — 1 a-nd 
for the spectra from panels (a), (c), (d) A is almost the same. 

4 Concluding remarks 

In this paper we apply the multifractal detrended fluctuation analysis (MFDFA) to the 
study of Drosophila ECG time series. On the example of Drosophila we have shown 



8 




Figure 4: The fractal spectrum h{q) for the time series of the parents, q changes contin- 
uously from —10 to 10. Panel (a): female parent. Panel (b): male parent 

that the presence of long-range correlations in the heartbeat activity is property not 
only of humans and complex animals and can be observed in much simpler animals as 
for example in Drosophila melanogaster. Opposite to the heartbeat dynamics of healtly 
humans which is described by broad range of Hurst exponents the intcrmaxima intervals 
of the time series of the heartbeat dynamics of healtly Drosophila do not have scaling 
properties and thus it cannot be described by means of scaling exponents and fractal 
spectra. We have shown that the presence of genetic defects can lead to long-range 
correlations of the heartbeat dynamics of Drosophila. The transfer of the multifractal 
properties from generation to generation and the similarity of the kinds and parameters 
of the multifractal spectra for different generations of Drosophila show that a correlation 
could exists between genetic properties and dynamic patterns in the heartbeat activity of 
simple animals like Drosophila. We can conjecture that the above correlation exists for 
the case of other simple animals and probably also for the case of more complex animals 
and ever humans. 
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Figure 5: The fractal spectrum h{q) for the time series of the first generation (the kids). 
From the top to the bottom the three characteristic shapes of this spectra are shown. 
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Figure 6: The fractal spectra h{q) for the time series of the second generation (the kids 
of the kids) . 
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parent; Panel (b): male parent. Parameters of the spectrum in panel (b) are: Q;^i„ = 0.35, 
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Figure 8: The fractal spectra f{a) for the time series of the second generation (the kids). 
Four characteristic shapes of this spectra are shown. The parameters of the spectra of 
parabohc kinds arc: Panel (a): «mm = 0.24, amax = 1-79, ai = 0.61, «r = 0.99, A = 0.28. 
fmax = 1-00 at a = 0.80. Panel (c): a^m = 0.56, amax = 0.72, ai = 0.56, ar = 0.65, 
A = 0.08 fmax = 1-00 at a = 0.59. Panel (d): a^in = 0.41, a^ax = 0.81, ai = 0.54, 
ar = 0.72, A = 0.18. f^ax = 1-00 at a = 0.65. 
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Figure 9: Parabolic fractal spectra f{ct) for the time series of the third generation (the 
kids of the kids). Four characteristic shapes of this spectra are shown. Parameters of 
the spectra are as follows: Panel (a): (Xmin — 0.44, amax — 0.71, ai — 0.55, ar — 0.67, 
A = 0.12. f^ax = 1.00 at a = 0.65. Panel (b): a™„ = 0.58, amax = 1.40, ai = 0.72, 
ar = 1.01, A = 0.29. f^ax = 1.00 at a = 0.86. Panel (c): a„i„ = 0.63, amax = 0.94, 
a = 0.72, ar = 0.83, A = 0.13./„„^ = 1.00 at a = 0.78. Panel (d): a^in = 0.12, 
(^max = 0.65, ai = 0.47, = 0.60, A = 0.13. fmax = LOO at a = 0.56. 
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